Biomolecular interactions determine how transcription factors recognize their DNA binding sites, how proteins interact with each other, and consequently how a biological system functions. Since many biological molecules bear considerable electric charge, electrostatic interactions are among the most important when studying biomolecular interactions. However, electrostatic interactions in biological systems are difficult to calculate accurately in practice. Aside from the significant charges carried by biomolecules such as DNA and proteins, the solvent itself, namely, water produces considerable electrostatic effects. Furthermore, hydrogen bonds, known to be involved in helix formation in both DNA and proteins, are essentially electrostatic in origin. Indeed, it seems that electrostatic effects often drive the physical-chemical processes in biological systems and, thereby, determine biological functions. Therefore, any attempt to perform molecular dynamics (MD) simulations of biological systems will require an adequate description of these electrostatic forces. Previously, we have developed a rigorous method to calculate the crucial electrostatic forces in a biomolecular system. Exact analytical results have been obtained for systems with sufficient symmetry. For example we have completely analyzed a system consisting of a pair of objects with planar surfaces containing embedded charges, which is an idealized model of molecular contact. We have also developed methods to compute the electrostatic forces for a biomolecular system in which the atoms are represented by spheres. The method is rigorous in the context of the model and the accuracy can be tuned to any desired level. However, the underlying idea we developed to compute the electrostatic forces in a biomolecular system can, with some modification, be applied to a more realistic and flexible model in which a biomolecule is represented by an arbitrary surface which is decomposed into a set of small triangular patches. There is, of course, fundamental effect that require more than electrostatics to understand. For example, the short distance electronic interactions call for quantum mechanical understanding. Recently, we have developed a quantum mechanical approach that is suitable for many electron systems by finding the density functional needed for the density functional theory. To see how we can integrate our new approach with existing ones, we have also try to use currently available density functionals to perform calculations for peptide fragmentation probability. The experience seems encouraging and some of the computational results were written and published. For the past year, we have been working intensively to test the accuracy and adequacy of surface triagularization by comparing the exact analytical result with those numerically obtained via surface triangularization. We found that the usual triangularization method (using flat triangles) is generally pathological. This pathology is somewhat expected as we know from classical electrostatics that pointy vertices and sharp edges will acquire infinite electric field. This pathological situation, plaguing all related numerical methods, disappears only if all the edges across all ajoining pathces are smooth. Using curved pathes such as spherical triangles and saddle surface pathes can in principle alleviate this problem. However, coming up with an exact analytic mathematical expression is very challenging. Nevertheless, we are determined to meet this challenge by finding the solution to this pathology.